Fit density models to cod and flounder an check residual co-occurence

Intro

In this script, I take the collated stomach data set and calculate aggregates (feeding ratio, total weight of prey groups) and predictor variables for diet data, aggregate to get 1 stomach = 1 row per prey type (not prey individual). I also select only the columns I need for model fitting, join environmental covariates and cpue covariates for cod and flounder, and lastly saduria biomass densities.

Load packages & source functions

# Load libraries, install if needed
library(tidyverse)
#> Warning: package 'tidyr' was built under R version 4.0.5
library(readxl)
library(tidylog)
library(RCurl)
library(RColorBrewer)
#> Warning: package 'RColorBrewer' was built under R version 4.0.5
library(patchwork)
library(janitor)
library(forcats)
library(gapminder)
library(viridis)
library(ggridges)
library(raster)
library(icesDatras)
library(ggalluvial)
library(ggrepel)
library(ncdf4)
library(chron)
library(rnaturalearth)
library(rnaturalearthdata)
library(mapplots)
library(geosphere)
library(quantreg)
#> Warning in .recacheSubclasses(def@className, def, env): undefined subclass
#> "numericVector" of class "Mnumeric"; definition not updated
library(brms)
#> Warning: package 'Rcpp' was built under R version 4.0.5
library(sdmTMB)
library(sp)
library(raster)
library(ggcorrplot)
library(ggpubr)
options(mc.cores = parallel::detectCores()) 

# Source code for map plots
source("/Users/maxlindmark/Dropbox/Max work/R/cod_interactions/R/functions/map_plot.R")

theme_set(theme_plot())

# Continuous colors
options(ggplot2.continuous.colour = "viridis")

# Discrete colors
scale_colour_discrete <- function(...) {
  scale_colour_brewer(palette = "Dark2")
}

scale_fill_discrete <- function(...) {
  scale_fill_brewer(palette = "Dark2")
}

Read and summarise data

cpue_full <- readr::read_csv("data/clean/catch_by_length_q1_q4.csv") %>%
  janitor::clean_names() %>% 
  rename(X = x,
         Y = y) 
#> Rows: 1759609 Columns: 22
#> ── Column specification ────────────────────────────────────────────────────────
#> Delimiter: ","
#> chr  (8): Country, haul.id, IDx, ices_rect, id_haul_stomach, species, haul.i...
#> dbl (14): density, year, lat, lon, quarter, Month, sub_div, length_cm, depth...
#> 
#> ℹ Use `spec()` to retrieve the full column specification for this data.
#> ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
#> rename: renamed 2 variables (X, Y)

# Summaries density by species group an haul
cpue <- cpue_full %>%
  mutate(species2 = species,
         species2 = ifelse(species == "cod" & length_cm >= 25, "large_cod", species),
         species2 = ifelse(species == "cod" & length_cm < 25, "small_cod", species2)) %>% 
  group_by(haul_id, species2) %>% 
  summarise(density = sum(density)) %>% 
  ungroup()
#> summarise: now 27,897 rows and 3 columns, one group variable remaining (haul_id)
#> ungroup: no grouping variables

# Trim data with quantiles
ggplot(cpue, aes(density)) + 
  geom_histogram() +
  facet_wrap(~species2, scales = "free")
#> `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

cpue <- cpue %>% 
  group_by(species2) %>% 
  mutate(dens_quants = quantile(density, probs = 0.995)) %>% 
  filter(density < dens_quants) %>% 
  ungroup() %>% 
  dplyr::select(-dens_quants)
#> filter (grouped): removed 140 rows (1%), 27,757 rows remaining
#> ungroup: no grouping variables

# Make wide
wcpue <- cpue %>% pivot_wider(names_from = species2, values_from = density)
#> pivot_wider: reorganized (species2, density) into (flounder, large_cod, small_cod) [was 27757x3, now 9372x4]

head(wcpue)
#> # A tibble: 6 × 4
#>   haul_id                  flounder large_cod small_cod
#>   <chr>                       <dbl>     <dbl>     <dbl>
#> 1 1993:1:GFR:SOL:H20:21:1     9.41      77.9    0.628  
#> 2 1993:1:GFR:SOL:H20:22:32    6.63       5.13   0      
#> 3 1993:1:GFR:SOL:H20:23:31    0.953      2.61   0.00459
#> 4 1993:1:GFR:SOL:H20:24:30    1.89       9.71   0      
#> 5 1993:1:GFR:SOL:H20:25:2    16.7      400.     4.78   
#> 6 1993:1:GFR:SOL:H20:26:3    13.5      179.     2.95

# Left join in trawl information
nrow(wcpue)
#> [1] 9372
hauls <- cpue_full %>% distinct(haul_id, .keep_all = TRUE) %>% dplyr::select(-density)
#> distinct: removed 1,750,236 rows (99%), 9,373 rows remaining
nrow(hauls)
#> [1] 9373

cpue2 <- left_join(wcpue, hauls)
#> Joining, by = "haul_id"
#> left_join: added 20 columns (year, lat, lon, quarter, country, …)
#> > rows only in x 0
#> > rows only in y ( 1)
#> > matched rows 9,372
#> > =======
#> > rows total 9,372

colnames(cpue2)
#>  [1] "haul_id"         "flounder"        "large_cod"       "small_cod"      
#>  [5] "year"            "lat"             "lon"             "quarter"        
#>  [9] "country"         "month"           "i_dx"            "ices_rect"      
#> [13] "sub_div"         "length_cm"       "id_haul_stomach" "species"        
#> [17] "haul_id_size"    "substrate"       "depth"           "temp"           
#> [21] "oxy"             "sal"             "X"               "Y"

cpue2 <- cpue2 %>% drop_na(oxy, temp, depth, sal, substrate, flounder, small_cod, large_cod)
#> drop_na: removed 506 rows (5%), 8,866 rows remaining

# Scale variables
cpue2 <- cpue2 %>% 
  mutate(depth_ct = depth - mean(depth),
         depth_sq = depth_ct * depth_ct,
         depth_sq_sc = (depth_sq - mean(depth_sq)) / sd(depth_sq),
         depth_sc = (depth - mean(depth)) / sd(depth),
         temp_ct = temp - mean(temp),
         temp_sq = temp_ct * temp_ct,
         temp_sq_sc = (temp_sq - mean(temp_sq)) / sd(temp_sq),
         temp_sc = (temp - mean(temp)) / sd(temp),
         oxy_sc = (oxy - mean(oxy)) / sd(oxy),
         sal_sc = (sal - mean(sal)) / sd(sal),
         fyear = as.factor(year),
         fsubstrate = as.factor(substrate))

Read and scale the prediction grid

pred_grid_1 <- read_csv("data/clean/pred_grid_(1_2).csv")
#> Rows: 296520 Columns: 18
#> ── Column specification ────────────────────────────────────────────────────────
#> Delimiter: ","
#> chr  (2): substrate, ices_rect
#> dbl (16): X, Y, depth, year, lon, lat, quarter, oxy, temp, sal, sub_div, den...
#> 
#> ℹ Use `spec()` to retrieve the full column specification for this data.
#> ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
pred_grid_2 <- read_csv("data/clean/pred_grid_(2_2).csv")
#> Rows: 296520 Columns: 18
#> ── Column specification ────────────────────────────────────────────────────────
#> Delimiter: ","
#> chr  (2): substrate, ices_rect
#> dbl (16): X, Y, depth, year, lon, lat, quarter, oxy, temp, sal, sub_div, den...
#> 
#> ℹ Use `spec()` to retrieve the full column specification for this data.
#> ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

pred_grid <- bind_rows(pred_grid_1, pred_grid_2)

# Scale with respect to data!
pred_grid <- pred_grid %>% 
  drop_na(substrate) %>% 
  mutate(X = X / 1000,
         Y = Y / 1000,
         depth_ct = depth - mean(cpue2$depth),
         depth_sq = depth_ct * depth_ct,
         depth_sq_sc = (depth_sq - mean(cpue2$depth_sq)) / sd(cpue2$depth_sq),
         depth_sc = (depth - mean(cpue2$depth)) / sd(cpue2$depth),
         temp_ct = temp - mean(cpue2$temp),
         temp_sq = temp_ct * temp_ct,
         temp_sq_sc = (temp_sq - mean(cpue2$temp_sq)) / sd(cpue2$temp_sq),
         temp_sc = (temp - mean(cpue2$temp)) / sd(cpue2$temp),
         oxy_sc = (oxy - mean(cpue2$oxy)) / sd(cpue2$oxy),
         sal_sc = (sal - mean(cpue2$sal)) / sd(cpue2$sal),
         fyear = as.factor(year),
         fsubstrate = as.factor(substrate))
#> drop_na: removed 280 rows (<1%), 592,760 rows remaining

# Split by quarter
pred_grid_q1 <- pred_grid %>% filter(quarter == 1)
#> filter: removed 296,380 rows (50%), 296,380 rows remaining
pred_grid_q4 <- pred_grid %>% filter(quarter == 4)
#> filter: removed 296,380 rows (50%), 296,380 rows remaining

# Load cache
# qwraps2::lazyload_cache_dir(path = "R/main_analysis/residual_cooccurence_cache/html")

Plot data in space

hist(cpue2$depth)

hist(log(cpue2$depth))

cpue_long <- cpue2 %>%
  pivot_longer(c(flounder, small_cod, large_cod)) %>% 
  rename(species2 = name, density = value) %>% 
  pivot_longer(c(depth, temp, oxy, sal)) %>% 
  rename(env_var = name, env_var_value = value)
#> pivot_longer: reorganized (flounder, large_cod, small_cod) into (name, value) [was 8866x36, now 26598x35]
#> rename: renamed 2 variables (species2, density)
#> pivot_longer: reorganized (depth, temp, oxy, sal) into (name, value) [was 26598x35, now 106392x33]
#> rename: renamed 2 variables (env_var, env_var_value)

ggplot(cpue_long, aes(density)) + 
  geom_histogram() + 
  facet_wrap(~species2, scales = "free", ncol = 1)
#> `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

ggplot(cpue_long, aes(env_var_value, density)) + 
  geom_point(alpha = 0.3) + 
  stat_smooth(method = "lm", formula = y ~ x, color = "blue", se = FALSE) + 
  stat_smooth(method = "lm", formula = y ~ x + I(x^2), color = "red", se = FALSE) + 
  facet_wrap(env_var~species2, scales = "free", ncol = 3) + 
  coord_cartesian(expand = 0)

cpue2 <- cpue2 %>% mutate(fle_presence = ifelse(flounder == 0, "N", "Y"),
                          l_cod_presence = ifelse(large_cod == 0, "N", "Y"),
                          s_cod_presence = ifelse(small_cod == 0, "N", "Y"))

# Biomass density
plot_map_fc + 
  geom_point(data = cpue2, aes(X*1000, Y*1000, color = flounder)) + 
  scale_color_viridis_c(trans = "sqrt") + 
  facet_wrap(~year)

plot_map_fc + 
  geom_point(data = cpue2, aes(X*1000, Y*1000, color = large_cod)) + 
  scale_color_viridis_c(trans = "sqrt") + 
  facet_wrap(~year)

plot_map_fc + 
  geom_point(data = cpue2, aes(X*1000, Y*1000, color = small_cod)) + 
  scale_color_viridis_c(trans = "sqrt") + 
  facet_wrap(~year)

# Presence / absence
plot_map_fc + 
  geom_point(data = cpue2, aes(X*1000, Y*1000, color = fle_presence)) + 
  facet_wrap(~year)

plot_map_fc + 
  geom_point(data = cpue2, aes(X*1000, Y*1000, color = l_cod_presence)) + 
  facet_wrap(~year)

plot_map_fc + 
  geom_point(data = cpue2, aes(X*1000, Y*1000, color = s_cod_presence)) + 
  facet_wrap(~year)

Create meshes

# Q1
spde_q1 <- make_mesh(filter(cpue2, quarter == 1),
                     xy_cols = c("X", "Y"),
                     n_knots = 150, 
                     type = "kmeans",
                     seed = 42)
#> filter: removed 3,595 rows (41%), 5,271 rows remaining

# All flounder Q4
spde_q4 <- make_mesh(filter(cpue2, quarter == 4),
                     xy_cols = c("X", "Y"),
                     n_knots = 150, 
                     type = "kmeans",
                     seed = 42)
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

Fit models

q1

# Small cod model q1
mcod_s_q1 <- sdmTMB(small_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                      temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                    data = filter(cpue2, quarter == 1),
                    mesh = spde_q1,
                    family = tweedie(link = "log"),
                    spatiotemporal = "off",
                    spatial = "off",
                    time = "year",
                    reml = FALSE,
                    control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 3,595 rows (41%), 5,271 rows remaining

sanity(mcod_s_q1)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
tidy(mcod_s_q1, conf.int = TRUE)
#>                             term    estimate  std.error   conf.low    conf.high
#> 1              fsubstratebedrock -1.96871848 0.91473061 -3.7615575 -0.175879429
#> 2            fsubstratehard clay  1.87652206 0.18105087  1.5216689  2.231375232
#> 3  fsubstratehard-bottom complex  0.24549231 0.24381343 -0.2323732  0.723357853
#> 4                  fsubstratemud  1.57164234 0.18474927  1.2095404  1.933744261
#> 5                 fsubstratesand  1.36482044 0.17999880  1.0120293  1.717611603
#> 6                      fyear1994 -0.89579165 0.24274470 -1.3715625 -0.420020777
#> 7                      fyear1995  0.95362905 0.21601896  0.5302397  1.377018437
#> 8                      fyear1996 -0.06513184 0.25416236 -0.5632809  0.433017233
#> 9                      fyear1997  0.22936794 0.22898119 -0.2194270  0.678162830
#> 10                     fyear1998  0.31542450 0.21378935 -0.1035949  0.734443921
#> 11                     fyear1999 -0.40932422 0.22221403 -0.8448557  0.026207275
#> 12                     fyear2000  0.60229507 0.22846295  0.1545159  1.050074225
#> 13                     fyear2001 -0.03191646 0.21118788 -0.4458371  0.382004187
#> 14                     fyear2002  1.34702437 0.21190442  0.9316993  1.762349396
#> 15                     fyear2003 -0.40537822 0.23154854 -0.8592050  0.048448587
#> 16                     fyear2004  0.14923236 0.20613502 -0.2547849  0.553249570
#> 17                     fyear2005  1.26278590 0.20097787  0.8688765  1.656695281
#> 18                     fyear2006  0.17871236 0.21764950 -0.2478728  0.605297537
#> 19                     fyear2007  0.48896819 0.20416217  0.0888177  0.889118688
#> 20                     fyear2008  0.81548756 0.20524194  0.4132208  1.217754371
#> 21                     fyear2009  0.45296670 0.20358098  0.0539553  0.851978091
#> 22                     fyear2010  1.22355963 0.20948396  0.8129786  1.634140645
#> 23                     fyear2011  0.10189229 0.21738472 -0.3241739  0.527958509
#> 24                     fyear2012  0.21738072 0.21079966 -0.1957790  0.630540467
#> 25                     fyear2013  1.01431599 0.20057926  0.6211879  1.407444112
#> 26                     fyear2014  0.78608232 0.20419583  0.3858658  1.186298800
#> 27                     fyear2015 -0.29174022 0.20535289 -0.6942245  0.110744053
#> 28                     fyear2016  0.25636094 0.20471823 -0.1448794  0.657601310
#> 29                     fyear2017  0.02954718 0.20554014 -0.3733041  0.432398459
#> 30                     fyear2018  0.21927612 0.20602953 -0.1845343  0.623086576
#> 31                     fyear2019  0.34945939 0.24231957 -0.1254782  0.824397024
#> 32                     fyear2020 -0.48783946 0.24593642 -0.9698660 -0.005812934
#> 33                      depth_sc  0.25653667 0.05349450  0.1516894  0.361383969
#> 34                   depth_sq_sc -0.97140657 0.03271048 -1.0355179 -0.907295212
#> 35                       temp_sc  0.40460077 0.07470586  0.2581800  0.551021562
#> 36                    temp_sq_sc -0.55014066 0.05308025 -0.6541760 -0.446105278
#> 37                        oxy_sc  0.66364798 0.05727951  0.5513822  0.775913752
#> 38                        sal_sc  0.49095133 0.03928212  0.4139598  0.567942870

# Small cod model q1 spatial
mcod_s_q1_space <- sdmTMB(small_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                            temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                          data = filter(cpue2, quarter == 1),
                          mesh = spde_q1,
                          family = tweedie(link = "log"),
                          spatiotemporal = "off",
                          spatial = "on",
                          time = "year",
                          reml = FALSE,
                          control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 3,595 rows (41%), 5,271 rows remaining

sanity(mcod_s_q1_space)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
#> ✔ Range parameter doesn't look unreasonably large
tidy(mcod_s_q1_space, conf.int = TRUE)
#>                             term    estimate  std.error     conf.low
#> 1              fsubstratebedrock -0.34285793 0.96268709 -2.229689962
#> 2            fsubstratehard clay  0.43350358 0.39738971 -0.345365938
#> 3  fsubstratehard-bottom complex  0.29473877 0.42795481 -0.544037245
#> 4                  fsubstratemud  0.39444265 0.39639261 -0.382472586
#> 5                 fsubstratesand  0.37238930 0.39674024 -0.405207280
#> 6                      fyear1994 -1.06671292 0.22761965 -1.512839234
#> 7                      fyear1995  0.82781201 0.20695939  0.422179053
#> 8                      fyear1996 -0.61654467 0.24232535 -1.091493628
#> 9                      fyear1997 -0.22146087 0.21674077 -0.646264979
#> 10                     fyear1998  0.16389757 0.20703306 -0.241879766
#> 11                     fyear1999 -0.54194991 0.21322485 -0.959862941
#> 12                     fyear2000  0.33309035 0.21463772 -0.087591844
#> 13                     fyear2001 -0.26585320 0.20211443 -0.661990216
#> 14                     fyear2002  1.15573938 0.20212193  0.759587682
#> 15                     fyear2003 -0.34426526 0.21369807 -0.763105787
#> 16                     fyear2004  0.14867751 0.19254607 -0.228705861
#> 17                     fyear2005  1.21009131 0.19384676  0.830158629
#> 18                     fyear2006 -0.24732020 0.20672130 -0.652486512
#> 19                     fyear2007  0.60918599 0.19582530  0.225375458
#> 20                     fyear2008  0.98563336 0.19706607  0.599390954
#> 21                     fyear2009  0.36128907 0.19775802 -0.026309530
#> 22                     fyear2010  0.92562453 0.20456613  0.524682285
#> 23                     fyear2011 -0.22851568 0.21804176 -0.655869686
#> 24                     fyear2012  0.05726625 0.20164832 -0.337957199
#> 25                     fyear2013  0.99936734 0.19435827  0.618432138
#> 26                     fyear2014  0.91231485 0.19478615  0.530541010
#> 27                     fyear2015 -0.02108612 0.18903748 -0.391592766
#> 28                     fyear2016  0.37927168 0.19244554  0.002085359
#> 29                     fyear2017  0.28487378 0.19388433 -0.095132527
#> 30                     fyear2018  0.28498202 0.20273812 -0.112377385
#> 31                     fyear2019  0.34139213 0.23511032 -0.119415629
#> 32                     fyear2020 -0.24123651 0.24569350 -0.722786930
#> 33                      depth_sc  1.07334553 0.09033358  0.896294969
#> 34                   depth_sq_sc -0.96795742 0.07084855 -1.106818033
#> 35                       temp_sc  0.23263628 0.08925010  0.057709296
#> 36                    temp_sq_sc -0.37949123 0.05343449 -0.484220910
#> 37                        oxy_sc  0.35652345 0.07343753  0.212588548
#> 38                        sal_sc  0.10284762 0.08248583 -0.058821649
#>      conf.high
#> 1   1.54397411
#> 2   1.21237309
#> 3   1.13351480
#> 4   1.17135789
#> 5   1.14998588
#> 6  -0.62058661
#> 7   1.23344496
#> 8  -0.14159571
#> 9   0.20334323
#> 10  0.56967491
#> 11 -0.12403687
#> 12  0.75377254
#> 13  0.13028381
#> 14  1.55189108
#> 15  0.07457526
#> 16  0.52606088
#> 17  1.59002398
#> 18  0.15784611
#> 19  0.99299653
#> 20  1.37187576
#> 21  0.74888767
#> 22  1.32656678
#> 23  0.19883832
#> 24  0.45248970
#> 25  1.38030254
#> 26  1.29408870
#> 27  0.34942052
#> 28  0.75645801
#> 29  0.66488009
#> 30  0.68234143
#> 31  0.80219989
#> 32  0.24031391
#> 33  1.25039608
#> 34 -0.82909681
#> 35  0.40756327
#> 36 -0.27476156
#> 37  0.50045836
#> 38  0.26451688

# Large cod model q1
mcod_l_q1 <- sdmTMB(large_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                      temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                    data = filter(cpue2, quarter == 1),
                    mesh = spde_q1,
                    family = tweedie(link = "log"),
                    spatiotemporal = "off",
                    spatial = "off",
                    time = "year",
                    reml = FALSE,
                    control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 3,595 rows (41%), 5,271 rows remaining

sanity(mcod_l_q1)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
tidy(mcod_l_q1, conf.int = TRUE)
#>                             term    estimate  std.error    conf.low   conf.high
#> 1              fsubstratebedrock  2.64514665 0.62883120  1.41266015  3.87763316
#> 2            fsubstratehard clay  4.95566879 0.13722582  4.68671113  5.22462646
#> 3  fsubstratehard-bottom complex  4.12395546 0.18142381  3.76837132  4.47953960
#> 4                  fsubstratemud  4.94273033 0.13873674  4.67081133  5.21464934
#> 5                 fsubstratesand  4.84618839 0.13597423  4.57968379  5.11269299
#> 6                      fyear1994  0.47773160 0.16922052  0.14606548  0.80939773
#> 7                      fyear1995  0.23560044 0.16895011 -0.09553569  0.56673657
#> 8                      fyear1996  1.43154880 0.17191402  1.09460351  1.76849409
#> 9                      fyear1997 -0.60124625 0.18216143 -0.95827609 -0.24421641
#> 10                     fyear1998 -0.62332593 0.17083282 -0.95815210 -0.28849975
#> 11                     fyear1999 -0.31706558 0.16811430 -0.64656356  0.01243240
#> 12                     fyear2000 -0.43787263 0.18499506 -0.80045629 -0.07528898
#> 13                     fyear2001 -0.09081981 0.16070454 -0.40579492  0.22415529
#> 14                     fyear2002  0.48599607 0.16697390  0.15873323  0.81325890
#> 15                     fyear2003  0.68291173 0.16611261  0.35733699  1.00848647
#> 16                     fyear2004 -0.56485537 0.16348225 -0.88527470 -0.24443604
#> 17                     fyear2005  0.12388484 0.15847307 -0.18671668  0.43448635
#> 18                     fyear2006  0.71466892 0.15934692  0.40235470  1.02698314
#> 19                     fyear2007  0.18574064 0.15558690 -0.11920409  0.49068537
#> 20                     fyear2008  0.60240283 0.15597332  0.29670074  0.90810492
#> 21                     fyear2009  0.91187997 0.15176386  0.61442827  1.20933167
#> 22                     fyear2010  1.41551882 0.15535173  1.11103503  1.72000262
#> 23                     fyear2011  1.60986737 0.15465254  1.30675395  1.91298078
#> 24                     fyear2012  0.85684664 0.15362888  0.55573956  1.15795371
#> 25                     fyear2013  0.76755378 0.15311557  0.46745278  1.06765479
#> 26                     fyear2014  0.03570462 0.15951161 -0.27693239  0.34834162
#> 27                     fyear2015 -0.06956322 0.15324438 -0.36991669  0.23079025
#> 28                     fyear2016  0.24319872 0.15343717 -0.05753260  0.54393004
#> 29                     fyear2017 -0.48692753 0.15692426 -0.79449343 -0.17936162
#> 30                     fyear2018  0.01985974 0.15586964 -0.28563916  0.32535863
#> 31                     fyear2019 -0.66374888 0.19117122 -1.03843758 -0.28906018
#> 32                     fyear2020 -0.84519285 0.18899363 -1.21561356 -0.47477214
#> 33                      depth_sc  0.18418004 0.03904185  0.10765942  0.26070066
#> 34                   depth_sq_sc -0.22902160 0.01684931 -0.26204564 -0.19599757
#> 35                       temp_sc  0.85122495 0.05433875  0.74472295  0.95772695
#> 36                    temp_sq_sc -0.27665011 0.03598891 -0.34718708 -0.20611314
#> 37                        oxy_sc  0.60631200 0.04426467  0.51955485  0.69306915
#> 38                        sal_sc  0.49208531 0.02893881  0.43536628  0.54880433

# Large cod model q1 spatial
mcod_l_q1_space <- sdmTMB(large_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                            temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                          data = filter(cpue2, quarter == 1),
                          mesh = spde_q1,
                          family = tweedie(link = "log"),
                          spatiotemporal = "off",
                          spatial = "on",
                          time = "year",
                          reml = FALSE,
                          control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 3,595 rows (41%), 5,271 rows remaining

sanity(mcod_l_q1_space)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
#> ✔ Range parameter doesn't look unreasonably large
tidy(mcod_l_q1_space, conf.int = TRUE)
#>                             term    estimate  std.error    conf.low   conf.high
#> 1              fsubstratebedrock  3.65588485 0.87486107  1.94118865  5.37058104
#> 2            fsubstratehard clay  4.05796820 0.61886634  2.84501246  5.27092394
#> 3  fsubstratehard-bottom complex  3.91807380 0.62896906  2.68531709  5.15083051
#> 4                  fsubstratemud  4.05533582 0.61860322  2.84289579  5.26777585
#> 5                 fsubstratesand  4.10477154 0.61857935  2.89237829  5.31716479
#> 6                      fyear1994  0.18125000 0.15387140 -0.12033241  0.48283240
#> 7                      fyear1995 -0.03392428 0.15435439 -0.33645333  0.26860477
#> 8                      fyear1996  0.68550527 0.16052079  0.37089030  1.00012025
#> 9                      fyear1997 -0.82301613 0.16705324 -1.15043447 -0.49559779
#> 10                     fyear1998 -0.95963374 0.15900091 -1.27126980 -0.64799768
#> 11                     fyear1999 -0.74958271 0.15826209 -1.05977071 -0.43939470
#> 12                     fyear2000 -0.81791543 0.17126527 -1.15358918 -0.48224167
#> 13                     fyear2001 -0.49120453 0.15019107 -0.78557363 -0.19683544
#> 14                     fyear2002  0.03790481 0.15543148 -0.26673530  0.34254492
#> 15                     fyear2003  0.20724132 0.15403867 -0.09466893  0.50915158
#> 16                     fyear2004 -0.91681063 0.15243842 -1.21558445 -0.61803681
#> 17                     fyear2005 -0.27054700 0.15004391 -0.56462766  0.02353365
#> 18                     fyear2006  0.27019300 0.14930842 -0.02244613  0.56283213
#> 19                     fyear2007 -0.06946656 0.14648117 -0.35656438  0.21763126
#> 20                     fyear2008  0.25007103 0.14828543 -0.04056307  0.54070513
#> 21                     fyear2009  0.53200891 0.14462146  0.24855606  0.81546177
#> 22                     fyear2010  0.73519650 0.14978083  0.44163146  1.02876153
#> 23                     fyear2011  0.68777729 0.15319983  0.38751114  0.98804344
#> 24                     fyear2012  0.31172507 0.14621413  0.02515064  0.59829951
#> 25                     fyear2013  0.15236074 0.14713985 -0.13602806  0.44074955
#> 26                     fyear2014 -0.24807200 0.14947917 -0.54104579  0.04490179
#> 27                     fyear2015  0.10082367 0.13978032 -0.17314072  0.37478805
#> 28                     fyear2016  0.12748695 0.14171746 -0.15027417  0.40524806
#> 29                     fyear2017 -0.39549093 0.14567352 -0.68100578 -0.10997607
#> 30                     fyear2018 -0.46970308 0.15293122 -0.76944277 -0.16996340
#> 31                     fyear2019 -0.76606776 0.17964960 -1.11817450 -0.41396103
#> 32                     fyear2020 -0.87080353 0.17973996 -1.22308738 -0.51851969
#> 33                      depth_sc  0.83137291 0.06439256  0.70516581  0.95758002
#> 34                   depth_sq_sc -0.30181568 0.02504008 -0.35089333 -0.25273803
#> 35                       temp_sc  0.46601213 0.06218866  0.34412459  0.58789968
#> 36                    temp_sq_sc -0.19308314 0.03535955 -0.26238660 -0.12377969
#> 37                        oxy_sc  0.22044214 0.05478142  0.11307253  0.32781175
#> 38                        sal_sc -0.19069159 0.05708485 -0.30257584 -0.07880734

# Flounder model q1
mfle_q1 <- sdmTMB(flounder ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                      temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                  data = filter(cpue2, quarter == 1),
                  mesh = spde_q1,
                  family = tweedie(link = "log"),
                  spatiotemporal = "off",
                  spatial = "off",
                  time = "year",
                  reml = FALSE,
                  control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 3,595 rows (41%), 5,271 rows remaining

sanity(mfle_q1)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
tidy(mfle_q1, conf.int = TRUE)
#>                             term    estimate  std.error     conf.low
#> 1              fsubstratebedrock  4.99666541 0.36333914  4.284533788
#> 2            fsubstratehard clay  4.05324009 0.13028270  3.797890689
#> 3  fsubstratehard-bottom complex  4.09193769 0.16132025  3.775755818
#> 4                  fsubstratemud  4.19828772 0.13122492  3.941091599
#> 5                 fsubstratesand  3.81330552 0.13144616  3.555675768
#> 6                      fyear1994 -0.24028243 0.15967160 -0.553233026
#> 7                      fyear1995 -0.05996836 0.15807077 -0.369781379
#> 8                      fyear1996 -0.09984072 0.16266080 -0.418650022
#> 9                      fyear1997 -0.87503465 0.16633438 -1.201044051
#> 10                     fyear1998 -0.74136915 0.15953873 -1.054059326
#> 11                     fyear1999 -0.45805245 0.15438894 -0.760649203
#> 12                     fyear2000  0.06366237 0.16614226 -0.261970478
#> 13                     fyear2001  0.16955433 0.14982675 -0.124100714
#> 14                     fyear2002  0.26903098 0.15729364 -0.039258892
#> 15                     fyear2003  0.27199168 0.15552675 -0.032835146
#> 16                     fyear2004 -0.72045920 0.15228548 -1.018933247
#> 17                     fyear2005 -0.10189060 0.14912847 -0.394177027
#> 18                     fyear2006  0.20256220 0.15070082 -0.092805976
#> 19                     fyear2007  0.29917698 0.14628949  0.012454848
#> 20                     fyear2008  0.34691978 0.14744652  0.057929913
#> 21                     fyear2009 -0.21485153 0.14704385 -0.503052174
#> 22                     fyear2010  0.22622153 0.14715017 -0.062187507
#> 23                     fyear2011 -0.13238407 0.14983484 -0.426054956
#> 24                     fyear2012 -0.14786599 0.14852436 -0.438968389
#> 25                     fyear2013 -0.05346312 0.14464140 -0.336955053
#> 26                     fyear2014  0.01602252 0.14938150 -0.276759832
#> 27                     fyear2015 -0.15947640 0.14547209 -0.444596456
#> 28                     fyear2016  0.08783420 0.14560139 -0.197539287
#> 29                     fyear2017 -0.12808688 0.14693267 -0.416069627
#> 30                     fyear2018  0.26403001 0.14645133 -0.023009319
#> 31                     fyear2019  0.15767834 0.16904385 -0.173641521
#> 32                     fyear2020  0.21928844 0.16875965 -0.111474389
#> 33                      depth_sc  0.52687056 0.03341193  0.461384388
#> 34                   depth_sq_sc -0.12804799 0.01281008 -0.153155289
#> 35                       temp_sc  0.31697398 0.04743742  0.223998341
#> 36                    temp_sq_sc  0.07181011 0.03301999  0.007092121
#> 37                        oxy_sc  0.33145736 0.03687451  0.259184643
#> 38                        sal_sc  0.11769582 0.02627909  0.066189739
#>      conf.high
#> 1   5.70879704
#> 2   4.30858949
#> 3   4.40811956
#> 4   4.45548384
#> 5   4.07093527
#> 6   0.07266816
#> 7   0.24984466
#> 8   0.21896859
#> 9  -0.54902525
#> 10 -0.42867898
#> 11 -0.15545569
#> 12  0.38929521
#> 13  0.46320937
#> 14  0.57732084
#> 15  0.57681851
#> 16 -0.42198515
#> 17  0.19039584
#> 18  0.49793037
#> 19  0.58589912
#> 20  0.63590964
#> 21  0.07334911
#> 22  0.51463057
#> 23  0.16128682
#> 24  0.14323642
#> 25  0.23002881
#> 26  0.30880487
#> 27  0.12564366
#> 28  0.37320769
#> 29  0.15989586
#> 30  0.55106935
#> 31  0.48899820
#> 32  0.55005127
#> 33  0.59235673
#> 34 -0.10294070
#> 35  0.40994962
#> 36  0.13652809
#> 37  0.40373008
#> 38  0.16920189

# Flounder model q1 spatial
mfle_q1_space <- sdmTMB(flounder ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                          temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                        data = filter(cpue2, quarter == 1),
                        mesh = spde_q1,
                        family = tweedie(link = "log"),
                        spatiotemporal = "off",
                        spatial = "on",
                        time = "year",
                        reml = FALSE,
                        control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 3,595 rows (41%), 5,271 rows remaining

sanity(mfle_q1_space)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
#> ✔ Range parameter doesn't look unreasonably large
tidy(mfle_q1_space, conf.int = TRUE)
#>                             term    estimate  std.error     conf.low
#> 1              fsubstratebedrock  4.43982047 0.44194014  3.573633712
#> 2            fsubstratehard clay  3.98176512 0.29154803  3.410341489
#> 3  fsubstratehard-bottom complex  3.95350012 0.30398390  3.357702623
#> 4                  fsubstratemud  4.08921272 0.29110860  3.518650354
#> 5                 fsubstratesand  3.92132335 0.29214697  3.348725818
#> 6                      fyear1994 -0.20426351 0.14365325 -0.485818711
#> 7                      fyear1995 -0.12674126 0.14279257 -0.406609545
#> 8                      fyear1996  0.14325487 0.15058892 -0.151893979
#> 9                      fyear1997 -0.90008026 0.15209145 -1.198174030
#> 10                     fyear1998 -0.88116454 0.14657044 -1.168437330
#> 11                     fyear1999 -0.43388223 0.14215160 -0.712494259
#> 12                     fyear2000 -0.07912842 0.15159357 -0.376246357
#> 13                     fyear2001  0.04078139 0.13828051 -0.230243424
#> 14                     fyear2002  0.27696207 0.14460591 -0.006460314
#> 15                     fyear2003  0.15635500 0.14125481 -0.120499347
#> 16                     fyear2004 -0.71543556 0.13849478 -0.986880347
#> 17                     fyear2005 -0.20875425 0.13706623 -0.477399130
#> 18                     fyear2006  0.21934440 0.13747581 -0.050103237
#> 19                     fyear2007  0.08649691 0.13446915 -0.177057779
#> 20                     fyear2008  0.28389764 0.13567208  0.017985252
#> 21                     fyear2009 -0.18510219 0.13579421 -0.451253949
#> 22                     fyear2010  0.48908851 0.13738793  0.219813116
#> 23                     fyear2011  0.17882479 0.14223733 -0.099955252
#> 24                     fyear2012 -0.09397562 0.13569689 -0.359936643
#> 25                     fyear2013  0.18001961 0.13374256 -0.082110990
#> 26                     fyear2014  0.08086310 0.13563375 -0.184974160
#> 27                     fyear2015 -0.31259885 0.13064713 -0.568662518
#> 28                     fyear2016  0.02571851 0.13098274 -0.231002933
#> 29                     fyear2017 -0.14272750 0.13307107 -0.403541996
#> 30                     fyear2018  0.22027872 0.13733562 -0.048894149
#> 31                     fyear2019  0.22194655 0.15656020 -0.084905807
#> 32                     fyear2020 -0.03263847 0.15585322 -0.338105171
#> 33                      depth_sc  0.19175149 0.05244555  0.088960090
#> 34                   depth_sq_sc -0.01460155 0.01933626 -0.052499914
#> 35                       temp_sc  0.57581452 0.05794221  0.462249874
#> 36                    temp_sq_sc -0.11362884 0.03270141 -0.177722432
#> 37                        oxy_sc  0.44660702 0.04601940  0.356410649
#> 38                        sal_sc  0.28956946 0.05497272  0.181824914
#>      conf.high
#> 1   5.30600722
#> 2   4.55318876
#> 3   4.54929762
#> 4   4.65977509
#> 5   4.49392089
#> 6   0.07729169
#> 7   0.15312703
#> 8   0.43840373
#> 9  -0.60198648
#> 10 -0.59389175
#> 11 -0.15527021
#> 12  0.21798952
#> 13  0.31180620
#> 14  0.56038445
#> 15  0.43320934
#> 16 -0.44399077
#> 17  0.05989063
#> 18  0.48879203
#> 19  0.35005160
#> 20  0.54981003
#> 21  0.08104958
#> 22  0.75836390
#> 23  0.45760482
#> 24  0.17198541
#> 25  0.44215022
#> 26  0.34670036
#> 27 -0.05653518
#> 28  0.28243996
#> 29  0.11808700
#> 30  0.48945160
#> 31  0.52879890
#> 32  0.27282823
#> 33  0.29454288
#> 34  0.02329682
#> 35  0.68937916
#> 36 -0.04953526
#> 37  0.53680339
#> 38  0.39731401

q4

# Small cod model q4
mcod_s_q4 <- sdmTMB(small_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                      temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                    data = filter(cpue2, quarter == 4),
                    mesh = spde_q4,
                    family = tweedie(link = "log"),
                    spatiotemporal = "off",
                    spatial = "off",
                    time = "year",
                    reml = FALSE,
                    control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

sanity(mcod_s_q4)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
tidy(mcod_s_q4, conf.int = TRUE)
#>                             term    estimate  std.error     conf.low
#> 1              fsubstratebedrock  0.24891861 0.52336315 -0.776854314
#> 2            fsubstratehard clay  0.92210091 0.22511990  0.480874021
#> 3  fsubstratehard-bottom complex -0.06514765 0.26950256 -0.593362960
#> 4                  fsubstratemud  0.44289773 0.22161580  0.008538744
#> 5                 fsubstratesand  0.49969302 0.22477919  0.059133898
#> 6                      fyear1994  1.29518505 0.27704358  0.752189603
#> 7                      fyear1995  0.41944485 0.30016066 -0.168859222
#> 8                      fyear1996  0.35633917 0.29028322 -0.212605491
#> 9                      fyear1997  0.96534678 0.27465703  0.427028895
#> 10                     fyear1998  0.51466119 0.28396588 -0.041901701
#> 11                     fyear1999  0.95559483 0.25754915  0.450807771
#> 12                     fyear2000  0.53111817 0.26929124  0.003317042
#> 13                     fyear2001  1.66771157 0.24652528  1.184530907
#> 14                     fyear2002  1.01481398 0.25228623  0.520342051
#> 15                     fyear2003  0.19424719 0.25754768 -0.310536983
#> 16                     fyear2004  1.81228339 0.24810180  1.326012793
#> 17                     fyear2005  1.03054074 0.24806946  0.544333534
#> 18                     fyear2006  1.24438646 0.24574263  0.762739741
#> 19                     fyear2007  1.22348055 0.24302929  0.747151902
#> 20                     fyear2008  0.67670204 0.24250016  0.201410464
#> 21                     fyear2009  1.32582989 0.23902614  0.857347264
#> 22                     fyear2010  0.98446906 0.24553347  0.503232310
#> 23                     fyear2011  0.21537810 0.24914744 -0.272941909
#> 24                     fyear2012  1.08422605 0.24989980  0.594431443
#> 25                     fyear2013  1.44332325 0.24177635  0.969450306
#> 26                     fyear2014  0.32076347 0.24887773 -0.167027928
#> 27                     fyear2015  0.34743992 0.24742583 -0.137505801
#> 28                     fyear2016 -0.15896761 0.24964373 -0.648260325
#> 29                     fyear2017  0.95798803 0.23638261  0.494686628
#> 30                     fyear2018 -0.45566045 0.25036403 -0.946364931
#> 31                     fyear2019 -0.24502349 0.28434075 -0.802321125
#> 32                     fyear2020 -0.40289823 0.28149186 -0.954612126
#> 33                      depth_sc  0.28619406 0.05327621  0.181774610
#> 34                   depth_sq_sc -1.30884731 0.04654316 -1.400070218
#> 35                       temp_sc  1.07958233 0.07529883  0.931999331
#> 36                    temp_sq_sc -0.28061097 0.04356847 -0.366003594
#> 37                        oxy_sc  0.65516641 0.06414230  0.529449814
#> 38                        sal_sc  0.13880384 0.04888545  0.042990116
#>      conf.high
#> 1   1.27469153
#> 2   1.36332781
#> 3   0.46306767
#> 4   0.87725672
#> 5   0.94025215
#> 6   1.83818049
#> 7   1.00774893
#> 8   0.92528384
#> 9   1.50366467
#> 10  1.07122409
#> 11  1.46038189
#> 12  1.05891929
#> 13  2.15089223
#> 14  1.50928590
#> 15  0.69903137
#> 16  2.29855399
#> 17  1.51674795
#> 18  1.72603317
#> 19  1.69980920
#> 20  1.15199362
#> 21  1.79431252
#> 22  1.46570582
#> 23  0.70369812
#> 24  1.57402065
#> 25  1.91719620
#> 26  0.80855486
#> 27  0.83238564
#> 28  0.33032510
#> 29  1.42128944
#> 30  0.03504403
#> 31  0.31227415
#> 32  0.14881568
#> 33  0.39061352
#> 34 -1.21762440
#> 35  1.22716533
#> 36 -0.19521834
#> 37  0.78088301
#> 38  0.23461757

# Small cod model q4 spatial
mcod_s_q4_space <- sdmTMB(small_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                            temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                          data = filter(cpue2, quarter == 4),
                          mesh = spde_q4,
                          family = tweedie(link = "log"),
                          spatiotemporal = "off",
                          spatial = "on",
                          time = "year",
                          reml = FALSE,
                          control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

sanity(mcod_s_q4_space)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✔ No fixed-effect standard errors look unreasonably large
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
#> ✔ Range parameter doesn't look unreasonably large
tidy(mcod_s_q4_space, conf.int = TRUE)
#>                             term    estimate  std.error     conf.low  conf.high
#> 1              fsubstratebedrock  1.07109590 0.64066171 -0.184577967  2.3267698
#> 2            fsubstratehard clay -0.20971994 0.33939805 -0.874927896  0.4554880
#> 3  fsubstratehard-bottom complex -0.41308682 0.36715815 -1.132703574  0.3065299
#> 4                  fsubstratemud -0.53657998 0.33706931 -1.197223699  0.1240637
#> 5                 fsubstratesand -0.23392951 0.33633748 -0.893138852  0.4252798
#> 6                      fyear1994  1.42137697 0.26766484  0.896763518  1.9459904
#> 7                      fyear1995  0.62955802 0.28495113  0.071064069  1.1880520
#> 8                      fyear1996  0.28985185 0.28395900 -0.266697567  0.8464013
#> 9                      fyear1997  1.20548438 0.26668264  0.682796013  1.7281727
#> 10                     fyear1998  0.64852108 0.27484400  0.109836739  1.1872054
#> 11                     fyear1999  1.15297787 0.25379155  0.655555571  1.6504002
#> 12                     fyear2000  0.77863883 0.26295138  0.263263591  1.2940141
#> 13                     fyear2001  2.12461229 0.24160333  1.651078453  2.5981461
#> 14                     fyear2002  1.43706789 0.24684350  0.953263518  1.9208723
#> 15                     fyear2003  0.74081014 0.25537432  0.240285677  1.2413346
#> 16                     fyear2004  2.36352114 0.24650853  1.880373292  2.8466690
#> 17                     fyear2005  1.45119795 0.24846229  0.964220801  1.9381751
#> 18                     fyear2006  1.72889893 0.24718511  1.244425024  2.2133728
#> 19                     fyear2007  1.85075126 0.24033986  1.379693791  2.3218087
#> 20                     fyear2008  1.26656927 0.24308519  0.790131058  1.7430075
#> 21                     fyear2009  1.77193727 0.23969689  1.302140005  2.2417345
#> 22                     fyear2010  1.44264919 0.24339115  0.965611296  1.9196871
#> 23                     fyear2011  0.67678496 0.24871753  0.189307571  1.1642624
#> 24                     fyear2012  1.55184638 0.24673385  1.068256916  2.0354358
#> 25                     fyear2013  1.88867237 0.24097675  1.416366624  2.3609781
#> 26                     fyear2014  0.84695989 0.25304037  0.351009892  1.3429099
#> 27                     fyear2015  0.91163475 0.25118034  0.419330321  1.4039392
#> 28                     fyear2016  0.43299230 0.25349364 -0.063846094  0.9298307
#> 29                     fyear2017  1.46681759 0.23847024  0.999424506  1.9342107
#> 30                     fyear2018  0.09662261 0.25387325 -0.400959817  0.5942050
#> 31                     fyear2019  0.26319964 0.28640090 -0.298135817  0.8245351
#> 32                     fyear2020  0.40694785 0.28219828 -0.146150624  0.9600463
#> 33                      depth_sc  0.18664716 0.09423981  0.001940531  0.3713538
#> 34                   depth_sq_sc -1.26823326 0.08632407 -1.437425340 -1.0990412
#> 35                       temp_sc  0.92072120 0.10677342  0.711449140  1.1299933
#> 36                    temp_sq_sc -0.29992427 0.05264981 -0.403116009 -0.1967325
#> 37                        oxy_sc  0.73336498 0.08310780  0.570476676  0.8962533
#> 38                        sal_sc  0.20983074 0.08908270  0.035231866  0.3844296

# Large cod model q4
mcod_l_q4 <- sdmTMB(large_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                      temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                    data = filter(cpue2, quarter == 4),
                    mesh = spde_q4,
                    family = tweedie(link = "log"),
                    spatiotemporal = "off",
                    spatial = "off",
                    time = "year",
                    reml = FALSE,
                    control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

sanity(mcod_l_q4)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✖ `b_j` standard error may be large
#> ℹ Try simplifying the model, adjusting the mesh, or adding priors
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
tidy(mcod_l_q4, conf.int = TRUE)
#>                             term     estimate  std.error     conf.low
#> 1              fsubstratebedrock  3.615410413 0.49180189  2.651496427
#> 2            fsubstratehard clay  4.955482953 0.16198501  4.637998170
#> 3  fsubstratehard-bottom complex  4.335607431 0.20623884  3.931386734
#> 4                  fsubstratemud  4.934428116 0.15697578  4.626761242
#> 5                 fsubstratesand  4.902833279 0.16571786  4.578032248
#> 6                      fyear1994  0.079547515 0.21371583 -0.339327808
#> 7                      fyear1995  0.497575775 0.21292310  0.080254173
#> 8                      fyear1996 -0.207014622 0.21674394 -0.631824929
#> 9                      fyear1997 -0.526577676 0.21649973 -0.950909345
#> 10                     fyear1998 -0.378526632 0.21451113 -0.798960726
#> 11                     fyear1999 -0.736796129 0.20244478 -1.133580611
#> 12                     fyear2000 -0.377177042 0.20334952 -0.775734786
#> 13                     fyear2001 -0.541367887 0.19371028 -0.921033056
#> 14                     fyear2002  0.289722959 0.18563492 -0.074114804
#> 15                     fyear2003 -1.048916073 0.19633190 -1.433719527
#> 16                     fyear2004 -0.601065281 0.19722507 -0.987619312
#> 17                     fyear2005  0.211219874 0.18218664 -0.145859386
#> 18                     fyear2006 -0.077320425 0.18378225 -0.437527016
#> 19                     fyear2007  0.160525298 0.18002683 -0.192320802
#> 20                     fyear2008  0.023261390 0.17858139 -0.326751704
#> 21                     fyear2009  0.002155868 0.17923441 -0.349137118
#> 22                     fyear2010  0.928157259 0.17468782  0.585775430
#> 23                     fyear2011  0.040940369 0.17798929 -0.307912222
#> 24                     fyear2012 -0.451508722 0.19098805 -0.825838418
#> 25                     fyear2013  0.023653636 0.18052082 -0.330160669
#> 26                     fyear2014 -0.087093528 0.18164254 -0.443106372
#> 27                     fyear2015  0.042492761 0.17807000 -0.306518017
#> 28                     fyear2016 -0.526261877 0.18348739 -0.885890562
#> 29                     fyear2017 -0.641298971 0.17665032 -0.987527234
#> 30                     fyear2018 -1.334124523 0.18908366 -1.704721685
#> 31                     fyear2019 -1.251554578 0.21766366 -1.678167508
#> 32                     fyear2020 -1.309280966 0.21314810 -1.727043572
#> 33                      depth_sc  0.547628725 0.04452566  0.460360026
#> 34                   depth_sq_sc -0.636524362 0.03525343 -0.705619813
#> 35                       temp_sc  0.831226263 0.06115285  0.711368877
#> 36                    temp_sq_sc -0.324120937 0.03502299 -0.392764729
#> 37                        oxy_sc  0.544387848 0.05217297  0.442130705
#> 38                        sal_sc  0.093280726 0.04318989  0.008630088
#>      conf.high
#> 1   4.57932440
#> 2   5.27296773
#> 3   4.73982813
#> 4   5.24209499
#> 5   5.22763431
#> 6   0.49842284
#> 7   0.91489738
#> 8   0.21779569
#> 9  -0.10224601
#> 10  0.04190746
#> 11 -0.34001165
#> 12  0.02138070
#> 13 -0.16170272
#> 14  0.65356072
#> 15 -0.66411262
#> 16 -0.21451125
#> 17  0.56829913
#> 18  0.28288617
#> 19  0.51337140
#> 20  0.37327449
#> 21  0.35344885
#> 22  1.27053909
#> 23  0.38979296
#> 24 -0.07717903
#> 25  0.37746794
#> 26  0.26891932
#> 27  0.39150354
#> 28 -0.16663319
#> 29 -0.29507071
#> 30 -0.96352736
#> 31 -0.82494165
#> 32 -0.89151836
#> 33  0.63489742
#> 34 -0.56742891
#> 35  0.95108365
#> 36 -0.25547715
#> 37  0.64664499
#> 38  0.17793136

# Large cod model q4 spatial
mcod_l_q4_space <- sdmTMB(large_cod ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                            temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                          data = filter(cpue2, quarter == 4),
                          mesh = spde_q4,
                          family = tweedie(link = "log"),
                          spatiotemporal = "off",
                          spatial = "on",
                          time = "year",
                          reml = FALSE,
                          control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

sanity(mcod_l_q4_space)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✖ `b_j` standard error may be large
#> ℹ Try simplifying the model, adjusting the mesh, or adding priors
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
#> ✔ Range parameter doesn't look unreasonably large
tidy(mcod_l_q4_space, conf.int = TRUE)
#>                             term     estimate  std.error    conf.low
#> 1              fsubstratebedrock  5.103730158 0.69328915  3.74490838
#> 2            fsubstratehard clay  4.589048584 0.49544560  3.61799305
#> 3  fsubstratehard-bottom complex  4.189598615 0.50983008  3.19035003
#> 4                  fsubstratemud  4.407270001 0.49451386  3.43804064
#> 5                 fsubstratesand  4.492862882 0.49609538  3.52053380
#> 6                      fyear1994 -0.099200696 0.18620052 -0.46414700
#> 7                      fyear1995  0.767425511 0.18087530  0.41291644
#> 8                      fyear1996 -0.381825002 0.18780855 -0.74992299
#> 9                      fyear1997 -0.433420044 0.18746204 -0.80083889
#> 10                     fyear1998 -0.331710939 0.18429710 -0.69292662
#> 11                     fyear1999 -0.826759473 0.17973058 -1.17902495
#> 12                     fyear2000 -0.276938741 0.17562043 -0.62114847
#> 13                     fyear2001 -0.481991023 0.16906600 -0.81335430
#> 14                     fyear2002  0.292109788 0.16247678 -0.02633885
#> 15                     fyear2003 -1.171889069 0.17754055 -1.51986215
#> 16                     fyear2004 -0.363447506 0.17431385 -0.70509638
#> 17                     fyear2005  0.079873466 0.16297020 -0.23954225
#> 18                     fyear2006 -0.214760846 0.16665803 -0.54140458
#> 19                     fyear2007  0.241032689 0.15986265 -0.07229235
#> 20                     fyear2008  0.252263694 0.16059012 -0.06248715
#> 21                     fyear2009  0.063806035 0.16023277 -0.25024442
#> 22                     fyear2010  0.533602252 0.15643821  0.22698900
#> 23                     fyear2011 -0.129273110 0.15829468 -0.43952499
#> 24                     fyear2012 -0.425513429 0.17047471 -0.75963772
#> 25                     fyear2013 -0.213562600 0.16164940 -0.53038959
#> 26                     fyear2014 -0.153749952 0.16743041 -0.48190752
#> 27                     fyear2015 -0.008369975 0.16300149 -0.32784702
#> 28                     fyear2016 -0.547449307 0.16781001 -0.87635089
#> 29                     fyear2017 -0.543151429 0.15747368 -0.85179418
#> 30                     fyear2018 -1.311951657 0.17301523 -1.65105528
#> 31                     fyear2019 -1.373505708 0.19900644 -1.76355117
#> 32                     fyear2020 -1.166575745 0.20079186 -1.56012056
#> 33                      depth_sc  0.430311101 0.07344864  0.28635441
#> 34                   depth_sq_sc -0.534917085 0.05787410 -0.64834824
#> 35                       temp_sc  0.077167176 0.07704927 -0.07384662
#> 36                    temp_sq_sc -0.015465077 0.03937601 -0.09264064
#> 37                        oxy_sc  0.715553722 0.06249348  0.59306875
#> 38                        sal_sc  0.306072548 0.06991218  0.16904720
#>      conf.high
#> 1   6.46255193
#> 2   5.56010412
#> 3   5.18884720
#> 4   5.37649936
#> 5   5.46519196
#> 6   0.26574561
#> 7   1.12193458
#> 8  -0.01372701
#> 9  -0.06600120
#> 10  0.02950474
#> 11 -0.47449400
#> 12  0.06727098
#> 13 -0.15062774
#> 14  0.61055843
#> 15 -0.82391599
#> 16 -0.02179864
#> 17  0.39928918
#> 18  0.11188288
#> 19  0.55435773
#> 20  0.56701454
#> 21  0.37785650
#> 22  0.84021550
#> 23  0.18097877
#> 24 -0.09138914
#> 25  0.10326439
#> 26  0.17440762
#> 27  0.31110707
#> 28 -0.21854772
#> 29 -0.23450868
#> 30 -0.97284804
#> 31 -0.98346025
#> 32 -0.77303093
#> 33  0.57426779
#> 34 -0.42148593
#> 35  0.22818097
#> 36  0.06171048
#> 37  0.83803870
#> 38  0.44309790

# Flounder model q4
mfle_q4 <- sdmTMB(flounder ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                      temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                  data = filter(cpue2, quarter == 4),
                  mesh = spde_q4,
                  family = tweedie(link = "log"),
                  spatiotemporal = "off",
                  spatial = "off",
                  time = "year",
                  reml = FALSE,
                  control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

sanity(mfle_q4)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✖ `b_j` standard error may be large
#> ℹ Try simplifying the model, adjusting the mesh, or adding priors
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
tidy(mfle_q4, conf.int = TRUE)
#>                             term      estimate  std.error    conf.low
#> 1              fsubstratebedrock  4.6615037843 0.31164341  4.05069392
#> 2            fsubstratehard clay  3.5854025940 0.16675474  3.25856931
#> 3  fsubstratehard-bottom complex  3.5480446664 0.19245757  3.17083476
#> 4                  fsubstratemud  3.4937256706 0.16173039  3.17673994
#> 5                 fsubstratesand  3.2696946138 0.16877245  2.93890668
#> 6                      fyear1994 -0.8045746409 0.23904040 -1.27308522
#> 7                      fyear1995 -0.0494297586 0.23062078 -0.50143818
#> 8                      fyear1996 -0.6111982199 0.23262630 -1.06713739
#> 9                      fyear1997 -0.4402297146 0.21813539 -0.86776722
#> 10                     fyear1998 -0.4878360037 0.22583882 -0.93047195
#> 11                     fyear1999 -0.9690916527 0.21803234 -1.39642718
#> 12                     fyear2000 -0.4733995937 0.21112683 -0.88720057
#> 13                     fyear2001 -0.7489172517 0.20733806 -1.15529238
#> 14                     fyear2002  0.0156969705 0.19249593 -0.36158812
#> 15                     fyear2003 -1.5798939628 0.20997303 -1.99143354
#> 16                     fyear2004 -0.4381240517 0.20453924 -0.83901359
#> 17                     fyear2005 -0.0576143994 0.18855510 -0.42717561
#> 18                     fyear2006 -0.3065833646 0.19088650 -0.68071403
#> 19                     fyear2007  0.0987308655 0.18641731 -0.26664034
#> 20                     fyear2008  0.2253096266 0.18166612 -0.13074943
#> 21                     fyear2009 -0.0241725482 0.18320105 -0.38324000
#> 22                     fyear2010  0.2437660843 0.18155547 -0.11207610
#> 23                     fyear2011  0.3711063583 0.18092545  0.01649899
#> 24                     fyear2012 -0.1890275344 0.19209616 -0.56552909
#> 25                     fyear2013  0.0910869113 0.18455673 -0.27063763
#> 26                     fyear2014 -0.4078437527 0.18874754 -0.77778214
#> 27                     fyear2015 -0.0251027775 0.18595258 -0.38956313
#> 28                     fyear2016 -0.2934659653 0.18372137 -0.65355323
#> 29                     fyear2017  0.1530211731 0.17964539 -0.19907732
#> 30                     fyear2018 -0.3269442525 0.18608949 -0.69167296
#> 31                     fyear2019 -0.3869720805 0.21084177 -0.80021435
#> 32                     fyear2020 -0.1989611616 0.20565603 -0.60203956
#> 33                      depth_sc -0.1720857706 0.04361497 -0.25756955
#> 34                   depth_sq_sc -0.5862771572 0.04628766 -0.67699930
#> 35                       temp_sc  0.0002634332 0.06090002 -0.11909842
#> 36                    temp_sq_sc  0.3807205063 0.03823918  0.30577309
#> 37                        oxy_sc  1.4861816951 0.05432021  1.37971603
#> 38                        sal_sc -0.0354741730 0.03983383 -0.11354704
#>      conf.high
#> 1   5.27231365
#> 2   3.91223588
#> 3   3.92525458
#> 4   3.81071141
#> 5   3.60048255
#> 6  -0.33606406
#> 7   0.40257867
#> 8  -0.15525905
#> 9  -0.01269221
#> 10 -0.04520006
#> 11 -0.54175613
#> 12 -0.05959861
#> 13 -0.34254212
#> 14  0.39298206
#> 15 -1.16835439
#> 16 -0.03723451
#> 17  0.31194681
#> 18  0.06754730
#> 19  0.46410207
#> 20  0.58136869
#> 21  0.33489490
#> 22  0.59960827
#> 23  0.72571373
#> 24  0.18747402
#> 25  0.45281145
#> 26 -0.03790537
#> 27  0.33935757
#> 28  0.06662129
#> 29  0.50511967
#> 30  0.03778445
#> 31  0.02627019
#> 32  0.20411724
#> 33 -0.08660199
#> 34 -0.49555501
#> 35  0.11962529
#> 36  0.45566792
#> 37  1.59264736
#> 38  0.04259869

# Flounder model q4 spatial
mfle_q4_space <- sdmTMB(flounder ~ 0 + fsubstrate + fyear + depth_sc + depth_sq_sc +
                          temp_sc + temp_sq_sc + oxy_sc + sal_sc,
                      data = filter(cpue2, quarter == 4),
                      mesh = spde_q4,
                      family = tweedie(link = "log"),
                      spatiotemporal = "off",
                      spatial = "on",
                      time = "year",
                      reml = FALSE,
                      control = sdmTMBcontrol(newton_loops = 1))
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

sanity(mfle_q4_space)
#> ✔ Non-linear minimizer suggests successful convergence
#> ✔ Hessian matrix is positive definite
#> ✔ No extreme or very small eigen values detected
#> ✔ No gradients with respect to fixed effects are >= 0.001
#> ✔ No fixed-effect standard errors are NA
#> ✖ `b_j` standard error may be large
#> ℹ Try simplifying the model, adjusting the mesh, or adding priors
#> ✔ No sigma parameters are < 0.01
#> ✔ No sigma parameters are > 100
#> ✔ Range parameter doesn't look unreasonably large
tidy(mfle_q4_space, conf.int = TRUE)
#>                             term    estimate  std.error     conf.low
#> 1              fsubstratebedrock  3.24794869 0.41960334  2.425541246
#> 2            fsubstratehard clay  2.86455220 0.33278839  2.212298941
#> 3  fsubstratehard-bottom complex  2.97553495 0.34265099  2.303951337
#> 4                  fsubstratemud  2.76045367 0.33122307  2.111268387
#> 5                 fsubstratesand  2.90546887 0.33207558  2.254612681
#> 6                      fyear1994 -0.79926186 0.20737851 -1.205716274
#> 7                      fyear1995  0.04604236 0.19743338 -0.340919947
#> 8                      fyear1996 -0.55170196 0.20322119 -0.950008165
#> 9                      fyear1997 -0.30983042 0.18976162 -0.681756361
#> 10                     fyear1998 -0.60026943 0.19699288 -0.986368392
#> 11                     fyear1999 -0.93999491 0.19717511 -1.326451035
#> 12                     fyear2000 -0.35205923 0.18783955 -0.720217973
#> 13                     fyear2001 -0.47545465 0.18423507 -0.836548742
#> 14                     fyear2002  0.42022487 0.16985311  0.087318890
#> 15                     fyear2003 -1.09229871 0.18788505 -1.460546635
#> 16                     fyear2004 -0.15233739 0.18284318 -0.510703436
#> 17                     fyear2005  0.11224188 0.17227607 -0.225413005
#> 18                     fyear2006 -0.01117985 0.17551077 -0.355174645
#> 19                     fyear2007  0.32808997 0.16823004 -0.001634851
#> 20                     fyear2008  0.40167141 0.16397614  0.080284086
#> 21                     fyear2009  0.21527348 0.16422729 -0.106606105
#> 22                     fyear2010  0.70075274 0.16200469  0.383229381
#> 23                     fyear2011  0.56831160 0.16077990  0.253188779
#> 24                     fyear2012  0.25415583 0.16850300 -0.076103976
#> 25                     fyear2013  0.39428580 0.16503076  0.070831462
#> 26                     fyear2014  0.16746106 0.17217498 -0.169995707
#> 27                     fyear2015  0.44613683 0.16857971  0.115726673
#> 28                     fyear2016  0.33316592 0.16601019  0.007791930
#> 29                     fyear2017  0.60390239 0.16055431  0.289221732
#> 30                     fyear2018  0.18781357 0.16992326 -0.145229908
#> 31                     fyear2019  0.18415294 0.19085396 -0.189913950
#> 32                     fyear2020  0.09858201 0.18470875 -0.263440495
#> 33                      depth_sc -0.13689884 0.06789870 -0.269977855
#> 34                   depth_sq_sc -0.62980393 0.06596533 -0.759093606
#> 35                       temp_sc  0.31635318 0.08319699  0.153290066
#> 36                    temp_sq_sc  0.07324120 0.04228867 -0.009643075
#> 37                        oxy_sc  1.28512872 0.06889794  1.150091252
#> 38                        sal_sc  0.20947782 0.07144799  0.069442340
#>       conf.high
#> 1   4.070356129
#> 2   3.516805454
#> 3   3.647118553
#> 4   3.409638961
#> 5   3.556325051
#> 6  -0.392807453
#> 7   0.433004666
#> 8  -0.153395748
#> 9   0.062095521
#> 10 -0.214170475
#> 11 -0.553538792
#> 12  0.016099521
#> 13 -0.114360554
#> 14  0.753130841
#> 15 -0.724050785
#> 16  0.206028664
#> 17  0.449896771
#> 18  0.332814939
#> 19  0.657814794
#> 20  0.723058737
#> 21  0.537153060
#> 22  1.018276090
#> 23  0.883434415
#> 24  0.584415634
#> 25  0.717740137
#> 26  0.504917822
#> 27  0.776546996
#> 28  0.658539900
#> 29  0.918583055
#> 30  0.520857047
#> 31  0.558219839
#> 32  0.460604512
#> 33 -0.003819832
#> 34 -0.500514251
#> 35  0.479416292
#> 36  0.156125470
#> 37  1.420166196
#> 38  0.349513305

Residuals in space

cpue2_q1 <- filter(cpue2, quarter == 1)
#> filter: removed 3,595 rows (41%), 5,271 rows remaining
cpue2_q4 <- filter(cpue2, quarter == 4)
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

cpue2_q1$res_small_cod <- residuals(mcod_s_q1)
cpue2_q1$res_large_cod <- residuals(mcod_l_q1)
cpue2_q1$res_flounder <- residuals(mfle_q1)

cpue2_q4$res_small_cod <- residuals(mcod_s_q4)
cpue2_q4$res_large_cod <- residuals(mcod_l_q4)
cpue2_q4$res_flounder <- residuals(mfle_q4)

cpue2_q1_long <- cpue2_q1 %>% pivot_longer(c(res_small_cod, res_large_cod, res_flounder))
#> pivot_longer: reorganized (res_small_cod, res_large_cod, res_flounder) into (name, value) [was 5271x42, now 15813x41]
cpue2_q4_long <- cpue2_q4 %>% pivot_longer(c(res_small_cod, res_large_cod, res_flounder))
#> pivot_longer: reorganized (res_small_cod, res_large_cod, res_flounder) into (name, value) [was 3595x42, now 10785x41]
  
res_long <- bind_rows(cpue2_q1_long, cpue2_q4_long) %>%
  rename(species2 = name) %>% 
  mutate(species2 = ifelse(species2 == "res_flounder", "Flounder", species2),
         species2 = ifelse(species2 == "res_large_cod", "Large cod", species2),
         species2 = ifelse(species2 == "res_small_cod", "Small cod", species2))
#> rename: renamed one variable (species2)

plot_map_fc + 
  geom_point(data = filter(res_long, year == 1999), aes(X*1000, Y*1000, color = value), size = 2) + 
  facet_grid(quarter ~ species2) + 
  scale_color_gradient2(name = "Residuals")
#> filter: removed 25,761 rows (97%), 837 rows remaining

ggsave("figures/supp/residuals_in_space.pdf", width = 20, height = 20, units = "cm")

Residual correlations

Raw correlation plot

# Plot correlation between variables
d_cor_q1 <- cpue2_q1 %>% dplyr::select(res_small_cod, res_large_cod, res_flounder)

d_cor_q4 <- cpue2_q4 %>% dplyr::select(res_small_cod, res_large_cod, res_flounder)

p1 <- ggplot(d_cor_q1, aes(res_large_cod, res_flounder)) + 
  geom_point(alpha = 0.2) + 
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red") +
  ggtitle("Quarter 1")

p2 <- ggplot(d_cor_q1, aes(res_small_cod, res_flounder)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Small cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red")

p3 <- ggplot(d_cor_q1, aes(res_large_cod, res_small_cod)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Small cod") +
  stat_cor(aes(label = ..r.label..), color = "red")

p4 <- ggplot(d_cor_q4, aes(res_large_cod, res_flounder)) + 
  geom_point(alpha = 0.2) + 
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red") +
  ggtitle("Quarter 4")

p5 <- ggplot(d_cor_q4, aes(res_small_cod, res_flounder)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Small cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red")

p6 <- ggplot(d_cor_q4, aes(res_large_cod, res_small_cod)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Small cod") +
  stat_cor(aes(label = ..r.label..), color = "red")

pp1 <- (p1 | p2 | p3) & theme(aspect.ratio = 1)
pp4 <- (p4 | p5 | p6) & theme(aspect.ratio = 1)

pp1 / pp4
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 2 rows containing non-finite values (stat_cor).
#> Warning: Removed 2 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).

ggsave("figures/raw_residuals.pdf", width = 22, height = 16, units = "cm")
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 2 rows containing non-finite values (stat_cor).
#> Warning: Removed 2 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).

corplot

corr_q1 <- round(cor(d_cor_q1), 3)
corr_q4 <- round(cor(d_cor_q4), 3)

pal <- brewer.pal(n = 3, name = "RdYlBu")

p1 <- ggcorrplot(corr_q1, hc.order = TRUE, type = "full", lab = TRUE, lab_size = 2.5,
                 #colors = c("steelblue", "white", "tomato")
                 colors = c(pal[3], "white", pal[1])) +
  theme_minimal(base_size = 11) + 
  labs(x = "", y = "") + 
  ggtitle("Quarter 1")

p4 <- ggcorrplot(corr_q4, hc.order = TRUE, type = "full", lab = TRUE, lab_size = 2.5,
                 #colors = c("steelblue", "white", "tomato")) +
                 colors = c(pal[3], "white", pal[1])) +
  theme_minimal(base_size = 11) + 
  labs(x = "", y = "") + 
  ggtitle("Quarter 4")

p1 + p4 + plot_layout(guides = "collect") & theme(legend.position = 'right')

ggsave("figures/supp/residu_corr.pdf", width = 20, height = 20, unit = "cm")

Residuals by subdiv

d_cor_q4 <- cpue2_q4 %>% dplyr::select(res_small_cod, res_large_cod, res_flounder, sub_div)

unique(d_cor_q4$sub_div)
#> [1] 24 25 27 28 26

d_cor_q4_24 <- d_cor_q4 %>% filter(sub_div == 24) %>% dplyr::select(-sub_div)
#> filter: removed 2,462 rows (68%), 1,133 rows remaining
d_cor_q4_25 <- d_cor_q4 %>% filter(sub_div == 25) %>% dplyr::select(-sub_div)
#> filter: removed 2,405 rows (67%), 1,190 rows remaining
d_cor_q4_26 <- d_cor_q4 %>% filter(sub_div == 26) %>% dplyr::select(-sub_div)
#> filter: removed 2,928 rows (81%), 667 rows remaining
d_cor_q4_27 <- d_cor_q4 %>% filter(sub_div == 27) %>% dplyr::select(-sub_div)
#> filter: removed 3,393 rows (94%), 202 rows remaining
d_cor_q4_28 <- d_cor_q4 %>% filter(sub_div == 28) %>% dplyr::select(-sub_div)
#> filter: removed 3,192 rows (89%), 403 rows remaining

corr_q4_24 <- round(cor(d_cor_q4_24), 3)
corr_q4_25 <- round(cor(d_cor_q4_25), 3)
corr_q4_26 <- round(cor(d_cor_q4_26), 3)
corr_q4_27 <- round(cor(d_cor_q4_27), 3)
corr_q4_28 <- round(cor(d_cor_q4_28), 3)

pal <- brewer.pal(n = 3, name = "RdYlBu")

p1 <- ggcorrplot(corr_q4_24, hc.order = TRUE, type = "full", lab = TRUE, lab_size = 2.5,
                 colors = c(pal[3], "white", pal[1])) +
  theme_minimal(base_size = 11) + 
  labs(x = "", y = "") + 
  ggtitle("Subdivision 24")
  
p2 <- ggcorrplot(corr_q4_25, hc.order = TRUE, type = "full", lab = TRUE, lab_size = 2.5,
                 colors = c(pal[3], "white", pal[1])) +
  theme_minimal(base_size = 11) + 
  labs(x = "", y = "") + 
  ggtitle("Subdivision 25")

p3 <- ggcorrplot(corr_q4_26, hc.order = TRUE, type = "full", lab = TRUE, lab_size = 2.5,
                 colors = c(pal[3], "white", pal[1])) +
  theme_minimal(base_size = 11) + 
  labs(x = "", y = "") + 
  ggtitle("Subdivision 26")


# p4 <- ggcorrplot(corr_q4_27, hc.order = TRUE, type = "full", lab = TRUE, lab_size = 2.5,
#                  colors = c(pal[3], "white", pal[1])) +
#   theme_minimal(base_size = 11) + 
#   labs(x = "", y = "") + 
#   ggtitle("Subdivision 27")

p5 <- ggcorrplot(corr_q4_28, hc.order = TRUE, type = "full", lab = TRUE, lab_size = 2.5,
                 colors = c(pal[3], "white", pal[1])) +
  theme_minimal(base_size = 11) + 
  labs(x = "", y = "") + 
  ggtitle("Subdivision 28")

#(p1 + p2 + p3) / (p4 + p5 + plot_spacer()) + plot_layout(guides = "collect") & theme(legend.position = 'bottom')
(p1 + p2 ) / (p3 + p5) + plot_layout(guides = "collect") & theme(legend.position = 'bottom')

ggsave("figures/supp/residu_corr_subdiv.pdf", width = 20, height = 20, unit = "cm")

cpue2 %>% 
  group_by(sub_div) %>% 
  summarise(n = n())
#> summarise: now 5 rows and 2 columns, ungrouped
#> # A tibble: 5 × 2
#>   sub_div     n
#>     <dbl> <int>
#> 1      24  2366
#> 2      25  3274
#> 3      26  1758
#> 4      27   431
#> 5      28  1037

raw correlation by ices rectangle!

d_cor_q4 <- cpue2_q4 %>%
  dplyr::select(res_small_cod,
                res_large_cod,
                res_flounder,
                #ices_rect,
                sub_div,
                flounder, 
                small_cod,
                large_cod,
                year) %>% 
  #mutate(year_rect = paste(year, ices_rect, sep = "_"))
  mutate(year_sd = paste(year, sub_div, sep = "_"))


d_cor_q4_sum <- cpue2_q4 %>%
  dplyr::select(res_small_cod,
                res_large_cod,
                res_flounder,
                #ices_rect,
                sub_div,
                flounder, 
                small_cod,
                large_cod,
                year) %>% 
  #mutate(year_rect = paste(year, ices_rect, sep = "_")) %>%
  mutate(year_sd = paste(year, sub_div, sep = "_")) %>% 
  #group_by(year_rect) %>%
  group_by(year_sd) %>% 
  summarise(res_small_cod = mean(res_small_cod),
            res_large_cod = mean(res_large_cod),
            res_flounder = mean(res_flounder),
            #ices_rect = head(ices_rect, 1),
            sub_div = head(sub_div, 1),
            flounder_dens = mean(log(flounder + 1)),
            small_cod_dens = mean(log(small_cod + 1)),
            large_cod_dens = mean(log(large_cod + 1)),
            year = head(year, 1)) %>% 
  ungroup() %>% 
  mutate(index = row_number())
#> summarise: now 138 rows and 9 columns, ungrouped
#> ungroup: no grouping variables

d_cor_q4_sum
#> # A tibble: 138 × 10
#>    year_sd res_sma…¹ res_l…² res_f…³ sub_div floun…⁴ small…⁵ large…⁶  year index
#>    <chr>       <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl> <dbl> <int>
#>  1 1993_24    0.145   0.158   0.443       24   3.28   1.02      4.48  1993     1
#>  2 1993_25   -0.718  -0.778  -1.41        25   0.665  0.440     2.04  1993     2
#>  3 1993_27    0.320   0.0350 -0.126       27   3.45   0.877     3.98  1993     3
#>  4 1993_28    0.0920 -0.810   0.0350      28   3.42   0.440     1.98  1993     4
#>  5 1994_24    0.170   0.0839  0.275       24   2.59   2.47      4.97  1994     5
#>  6 1994_25   -0.665   0.661  -0.0524      25   1.33   0.969     4.60  1994     6
#>  7 1994_26   -0.104   0.212  -0.750       26   0      0.0233    2.71  1994     7
#>  8 1994_27   -0.751  -1.08   -0.634       27   0.740  0         1.09  1994     8
#>  9 1994_28   -0.916  -0.689   0.843       28   3.39   0.181     2.47  1994     9
#> 10 1995_24    0.181   0.277   0.218       24   3.22   1.66      5.75  1995    10
#> # … with 128 more rows, and abbreviated variable names ¹​res_small_cod,
#> #   ²​res_large_cod, ³​res_flounder, ⁴​flounder_dens, ⁵​small_cod_dens,
#> #   ⁶​large_cod_dens

data_list <- list()

#for(i in unique(d_cor_q4_sum$year_rect)){
for(i in unique(d_cor_q4_sum$year_sd)){

  # dd <- filter(d_cor_q4, year_rect == i)
  # ddd <- filter(d_cor_q4_sum, year_rect == i)
  
  dd <- filter(d_cor_q4, year_sd == i)
  ddd <- filter(d_cor_q4_sum, year_sd == i)
  
  cor_df <- round(cor(dd[1:3]), 3)
  
  cor_df2 <- data.frame(small_cod_large_cod = cor_df[1, 2],
                        small_cod_flounder = cor_df[1, 3],
                        large_cod_flounder = cor_df[2, 3])
  
  #cor_df2$year_rect <- i
  cor_df2$year_sd <- i
  cor_df2$year <- ddd$year
  #cor_df2$ices_rect <- ddd$ices_rect
  cor_df2$sub_div <- ddd$sub_div
  cor_df2$flounder <- ddd$flounder_dens
  cor_df2$small_cod <- ddd$small_cod_dens
  cor_df2$large_cod <- ddd$large_cod_dens
  
  data_list[[ddd$index]] <- cor_df2

}
#> filter: removed 3,560 rows (99%), 35 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,582 rows (>99%), 13 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,588 rows (>99%), 7 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,558 rows (99%), 37 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,586 rows (>99%), 9 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,592 rows (>99%), 3 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,588 rows (>99%), 7 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,590 rows (>99%), 5 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,565 rows (99%), 30 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,586 rows (>99%), 9 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,592 rows (>99%), 3 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,588 rows (>99%), 7 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,558 rows (99%), 37 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,581 rows (>99%), 14 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,593 rows (>99%), 2 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,558 rows (99%), 37 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,578 rows (>99%), 17 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,589 rows (>99%), 6 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,584 rows (>99%), 11 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,561 rows (99%), 34 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,580 rows (>99%), 15 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,592 rows (>99%), 3 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,590 rows (>99%), 5 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,585 rows (>99%), 10 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,553 rows (99%), 42 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,557 rows (99%), 38 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,589 rows (>99%), 6 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,556 rows (99%), 39 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,570 rows (99%), 25 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,584 rows (>99%), 11 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,588 rows (>99%), 7 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,588 rows (>99%), 7 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,553 rows (99%), 42 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,545 rows (99%), 50 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,587 rows (>99%), 8 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,586 rows (>99%), 9 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,546 rows (99%), 49 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,551 rows (99%), 44 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,583 rows (>99%), 12 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,585 rows (>99%), 10 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,551 rows (99%), 44 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,552 rows (99%), 43 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,580 rows (>99%), 15 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,590 rows (>99%), 5 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,576 rows (99%), 19 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,550 rows (99%), 45 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,569 rows (99%), 26 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,579 rows (>99%), 16 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,592 rows (>99%), 3 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,581 rows (>99%), 14 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,546 rows (99%), 49 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,566 rows (99%), 29 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,557 rows (99%), 38 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,589 rows (>99%), 6 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,575 rows (99%), 20 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,550 rows (99%), 45 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,550 rows (99%), 45 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,557 rows (99%), 38 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,590 rows (>99%), 5 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,579 rows (>99%), 16 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,552 rows (99%), 43 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,545 rows (99%), 50 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,555 rows (99%), 40 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,587 rows (>99%), 8 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,575 rows (99%), 20 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,545 rows (99%), 50 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,553 rows (99%), 42 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,549 rows (99%), 46 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,585 rows (>99%), 10 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,570 rows (99%), 25 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,538 rows (98%), 57 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,537 rows (98%), 58 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,564 rows (99%), 31 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,588 rows (>99%), 7 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,578 rows (>99%), 17 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,565 rows (99%), 30 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,536 rows (98%), 59 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,555 rows (99%), 40 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,590 rows (>99%), 5 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,573 rows (99%), 22 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,541 rows (98%), 54 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,537 rows (98%), 58 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,562 rows (99%), 33 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,584 rows (>99%), 11 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,576 rows (99%), 19 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,549 rows (99%), 46 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,569 rows (99%), 26 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,564 rows (99%), 31 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,585 rows (>99%), 10 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,578 rows (>99%), 17 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,555 rows (99%), 40 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,533 rows (98%), 62 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,575 rows (99%), 20 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,590 rows (>99%), 5 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,575 rows (99%), 20 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,542 rows (99%), 53 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,531 rows (98%), 64 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,566 rows (99%), 29 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,586 rows (>99%), 9 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,575 rows (99%), 20 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,554 rows (99%), 41 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,529 rows (98%), 66 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,564 rows (99%), 31 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,586 rows (>99%), 9 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,578 rows (>99%), 17 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,540 rows (98%), 55 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,516 rows (98%), 79 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,543 rows (99%), 52 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,585 rows (>99%), 10 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,579 rows (>99%), 16 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,544 rows (99%), 51 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,499 rows (97%), 96 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,549 rows (99%), 46 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,587 rows (>99%), 8 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,573 rows (99%), 22 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,552 rows (99%), 43 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,497 rows (97%), 98 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,562 rows (99%), 33 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,586 rows (>99%), 9 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,568 rows (99%), 27 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,594 rows (>99%), one row remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,561 rows (99%), 34 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,555 rows (99%), 40 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,585 rows (>99%), 10 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,590 rows (>99%), 5 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,591 rows (>99%), 4 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,574 rows (99%), 21 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,559 rows (99%), 36 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,583 rows (>99%), 12 rows remaining
#> filter: removed 137 rows (99%), one row remaining
#> filter: removed 3,567 rows (99%), 28 rows remaining
#> filter: removed 137 rows (99%), one row remaining

big_dat <- dplyr::bind_rows(data_list)

big_dat$flounder_plus_small_cod <- big_dat$flounder + big_dat$small_cod
big_dat$flounder_plus_large_cod <- big_dat$flounder + big_dat$large_cod

ggplot(filter(big_dat, small_cod_flounder > -1 & small_cod_flounder < 1),
       aes(flounder_plus_small_cod, small_cod_flounder)) + 
  geom_point()
#> filter: removed 2 rows (1%), 136 rows remaining

ggplot(filter(big_dat, large_cod_flounder > -1 & large_cod_flounder < 1),
       aes(flounder_plus_large_cod, large_cod_flounder)) + 
  geom_point()
#> filter: removed 3 rows (2%), 135 rows remaining

mcod_s_q1_pred <- predict(mcod_s_q1_space, newdata = pred_grid_q1)
mcod_l_q1_pred <- predict(mcod_l_q1_space, newdata = pred_grid_q1)
mfle_q1_pred <- predict(mfle_q1_space, newdata = pred_grid_q1)

mcod_s_q4_pred <- predict(mcod_s_q4_space, newdata = pred_grid_q4)
mcod_l_q4_pred <- predict(mcod_l_q4_space, newdata = pred_grid_q4)
mfle_q4_pred <- predict(mfle_q4_space, newdata = pred_grid_q4)

omega_plot <- bind_rows(mcod_s_q1_pred %>% mutate(species = "Cod < 25 cm"),
                        mcod_s_q4_pred %>% mutate(species = "Cod < 25 cm"),
                        mcod_l_q1_pred %>% mutate(species = "Cod > 25 cm"),
                        mcod_l_q4_pred %>% mutate(species = "Cod > 25 cm"),
                        mcod_l_q1_pred %>% mutate(species = "Flounder"),
                        mcod_l_q4_pred %>% mutate(species = "Flounder"))


plot_map +
  geom_raster(data = filter(omega_plot, year == 1999),
              aes(X*1000, Y*1000, fill = omega_s)) +
  facet_grid(quarter ~ species) +
  scale_fill_gradient2()
#> filter: removed 1,714,770 rows (96%), 63,510 rows remaining

ggsave("figures/supp/omega_random_field.pdf", width = 20, height = 20, units = "cm")

# Now plot residual correlation
cpue2_spatial_q1 <- filter(cpue2, quarter == 1)
#> filter: removed 3,595 rows (41%), 5,271 rows remaining
cpue2_spatial_q4 <- filter(cpue2, quarter == 4)
#> filter: removed 5,271 rows (59%), 3,595 rows remaining

cpue2_spatial_q1$res_small_cod <- residuals(mcod_s_q1_space)
cpue2_spatial_q1$res_large_cod <- residuals(mcod_l_q1_space)
cpue2_spatial_q1$res_flounder <- residuals(mfle_q1_space)

cpue2_spatial_q4$res_small_cod <- residuals(mcod_s_q4_space)
cpue2_spatial_q4$res_large_cod <- residuals(mcod_l_q4_space)
cpue2_spatial_q4$res_flounder <- residuals(mfle_q4_space)

# Plot correlation between variables
d_cor_q1_spatial <- cpue2_spatial_q1 %>% dplyr::select(res_small_cod,
                                                       res_large_cod, 
                                                       res_flounder)

d_cor_q4_spatial <- cpue2_spatial_q1 %>% dplyr::select(res_small_cod,
                                                       res_large_cod,
                                                       res_flounder)

# min(d_cor_q4$res_small_cod, d_cor_q4$res_large_cod, d_cor_q4$res_flounder)
# max(d_cor_q4$res_small_cod, d_cor_q4$res_large_cod, d_cor_q4$res_flounder)

p1 <- ggplot(d_cor_q1_spatial, aes(res_large_cod, res_flounder)) + 
  geom_point(alpha = 0.2) + 
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red") +
  ggtitle("Quarter 1")

p2 <- ggplot(d_cor_q1_spatial, aes(res_small_cod, res_flounder)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Small cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red")

p3 <- ggplot(d_cor_q1_spatial, aes(res_large_cod, res_small_cod)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Small cod") +
  stat_cor(aes(label = ..r.label..), color = "red")

p4 <- ggplot(d_cor_q4_spatial, aes(res_large_cod, res_flounder)) + 
  geom_point(alpha = 0.2) + 
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red") + 
  ggtitle("Quarter 4")

p5 <- ggplot(d_cor_q4_spatial, aes(res_small_cod, res_flounder)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Small cod", y = "Flounder") +
  stat_cor(aes(label = ..r.label..), color = "red")

p6 <- ggplot(d_cor_q4_spatial, aes(res_large_cod, res_small_cod)) + 
  geom_point(alpha = 0.2) +
  xlim(-6, 6) +
  ylim(-6, 6) +
  labs(x = "Large cod", y = "Small cod") +
  stat_cor(aes(label = ..r.label..), color = "red")

pp1 <- (p1 | p2 | p3) & theme(aspect.ratio = 1)
pp4 <- (p4 | p5 | p6) & theme(aspect.ratio = 1)

(pp1 / pp4)
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).

ggsave("figures/supp/raw_residuals_spatial_model.pdf", width = 22, height = 16, units = "cm")
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).
#> Warning: Removed 1 rows containing non-finite values (stat_cor).
#> Warning: Removed 1 rows containing missing values (geom_point).